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ABSTRACT 

The  authors  investigate  the  effect  of  a  change  in  the  rate  of  thermal  damping  upon  the  climate  of  an  ocean 
general  circulation  model.  Initially,  the  thermal  forcing  condition  is  that  proposed  by  Haney,  that  is,  restoring 
the  model  surface  temperature  to  a  climatology.  The  restoring  condition  represents  a  strong  damping.  When  a 
steady  state  is  reached,  the  thermal  damping  is  switched  to  a  weaker  one,  but  the  atmosphere -ocean  heat 
exchange  is  adjusted  so  that  at  the  moment  of  the  switch  the  heat  flux  is  identical  to  that  prior  to  the  switch.  It 
is  found  that  interdecadal  oscillations  and  climate  drift  occur  as  a  result  of  the  switch,  regardless  of  the  forcing 
condition  for  salinity.  The  cause  for  the  variability  and  drift  can  be  traced  to  the  spinup.  During  the  spinup,  the 
surface  climatology  of  the  model  ocean  is  forcefully  “nudged”  toward  that  of  the  climatology,  regardless  of 
whether  or  not  the  internal  dynamics  of  the  model  ocean  can  maintain  the  climatology.  This  leads  to  intermittent 
convections  in  the  spinup  state.  When  the  thermal  damping  becomes  weaker,  the  system  chooses  a  convective 
pattern  (the  location  and  intensity  of  the  convection)  more  compatible  with  the  internal  dynamics.  An  implication 
of  these  results  is  that  drift  and  variability  in  a  coupled  model  may  be  caused  by  the  mechanism.  Effects  of  flux 
corrections  in  coupled  models  are  discussed. 


1.  Introduction 

In  atmosphere-ocean  coupling,  a  precoupling 
spinup  is  often  carried  out  to  provide  initial  conditions 
and  to  derive  flux  corrections.  During  this  spinup, 
strong  restoring  forces  on  both  temperature  and  salinity 
are  usually  used  in  order  to  best  represent  the  current 
ocean  climate.  These  represent  strong  thermal  and  sa¬ 
line  damping.  (Ocean  modelers  often  refer  to  such  a 
procedure  as  a  strong  coupling  because  the  surface 
ocean  conditions  are  strongly  constrained  by  the  at¬ 
mospheric  conditions  represented  by  the  restoring  forc¬ 
ings.)  Upon  switching  from  restoring  boundary  condi¬ 
tions  to  an  active  atmosphere,  the  damping  through  sa¬ 
line  forcing  weakens  and  the  interaction  between 
freshwater  flux  and  sea  surface  salinity  (SSS)  becomes 
indirect  and  nonlinear,  mostly  via  the  dependence  of 
evaporation  upon  sea  surface  temperature. 

Based  on  such  a  change,  many  studies  have  ex¬ 
amined  the  effect  of  a  change  of  the  saline  forcing 
from  a  restoring  to  a  flux  condition,  while  maintain¬ 
ing  the  restoration  for  temperature  (e.g.,  Bryan  1986; 
Marotzke  1990;  Marotzke  and  Willebrand  1991; 
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Weaver  and  Sarachik  1991a,b;  Moore  and  Reason 
1993;  Zhang  et  al.  1993;  Power  and  Kleeman  1993; 
Power  et  al.  1994;  Huang  1994;  Power  1995;  Cai  and 
Godfrey  1995;  Cai  1995a,b).  This  set  of  forcing  con¬ 
ditions  is  often  called  “mixed  boundary  conditions.” 
These  previous  studies  have  revealed  interesting 
features  of  the  thermohaline  circulation,  including  a 
polar  halocline  catastrophe  (Bryan  1986),  a  flush 
(Marotzke  1989;  Weaver  and  Sarachik  1991b),  low- 
frequency  oscillations  (Marotzke  1989),  and  oscil¬ 
lations  on  decadal  timescale  (Weaver  and  Sarachik 
1991a,b;  Cai  1995a).  Further,  the  thermohaline  cir¬ 
culation  is  extremely  sensitive  to  perturbations  in  the 
freshwater  flux.  A  small  perturbation  can  “flip”  the 
ocean  circulation  from  an  equilibrium  with  a  north¬ 
ern  sinking  cell  to  one  without  (Marotzke  1990;  Mar¬ 
otzke  and  Willebrand  1991;  Hughes  and  Weaver 
1994),  and  a  substantial  drift  can  be  induced  by  noise 
in  the  freshwater  flux  (Power  1995). 

It  is  now  clear  that  when  the  thermal  forcing  is 
switched  to  that  provided  by  an  active  atmosphere  the 
thermal  damping  is  also  weakened  (Zhang  et  al.  1993; 
Moore  and  Reason  1993;  Mikolajewicz  and  Maier-Rei- 
mer  1994;  Marotzke  1994;  Rahmstorf  and  Willebrand 
1995;  Cai  and  Godfrey  1995).  In  particular,  Marotzke 
(1994)  suggested  that  the  thermal  damping  rate  is 
about  12  W  m-2  K-1,  rather  than  the  commonly  used 
30-50  W  m-2  K'1.  This  is  confirmed  by  Rahmstorf 
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and  Willebrand  (1995)  in  an  ocean-only  model  forced 
by  a  thermal-energy-balanced  atmosphere  model. 
Hence,  it  is  equally  important  to  study  the  effect  of  a 
change  in  thermal  damping  because  all  coupled  mod¬ 
els,  flux  corrected  or  not,  go  through  such  a  change 
when  coupling  occurs. 

Most  of  the  coupled  atmosphere-ocean-sea  ice 
models  (e.g.,  Manabe  and  Stouffer  1993,  1994;  Cu- 
basch  et  al.  1994;  Gordon  and  O’ Farrell  1996,  hereafter 
GOF)  employ  a  “flux  correction”  technique  (Sausen 
et  al.  1988;  Sausen  and  Lunkeit  1990;  Moore  and  Gor¬ 
don  1994)  to  avoid  “climate  drift,”  which  is  the  trend 
of  each  model  component  with  time  after  coupling,  to¬ 
ward  a  climatology  significantly  different  from  that 
prior  to  coupling.  The  flux  correction  fields  to  the  ocean 
model  are  usually  obtained  and  applied  in  the  following 
way.  The  ocean  component  is  spun  up  under  restoring 
boundary  conditions  (Haney  1971),  which  force  the 
model  toward  SST  and  SSS  climatologies,  and  the  im¬ 
plied  fluxes  required  for  maintaining  these  climatolo¬ 
gies  are  diagnosed.  Similarly,  the  atmospheric  com¬ 
ponent  is  spun  up  under  the  influence  of  the  same  cli¬ 
matologies,  and  another  set  of  implied  freshwater  and 
heat  exchange  with  the  ocean  is  diagnosed.  The  differ¬ 
ence  between  these  fluxes  is  the  flux  correction;  this  is 
applied  during  subsequent  coupled  experiments.  This 
technique  in  effect  ensures  that  at  the  moment  of  cou¬ 
pling  the  ocean  receives  heat  and  freshwater  fluxes  ap¬ 
proximately  the  same  as  those  prior  to  the  coupling. 
Thereafter,  the  total  flux  is  varying  according  to  the 
dynamics  of  the  coupled  system. 

However,  the  correction  does  not  fully  prevent  drift, 
particularly  at  high  latitudes  (e.g.,  Figs.  5  and  6  of 
GOF).  Further,  upon  coupling,  oceanic  variabilities  of 
interdecadal  timescale  occur  (Delworth  et  al.  1993). 
An  examination  of  this  issue  requires  an  understanding 
as  to  what  may  induce  such  oscillations  and  anomalous 
drift.  In  this  study,  we  suggest  that  a  possible  cause  for 
drift  and  variability  is  the  change  in  the  rate  of  thermal 
damping. 

The  structure  for  the  rest  of  this  paper  is  as  follows. 
In  section  2,  the  ocean  general  circulation  model 
(OGCM)  and  the  thermal  forcing  condition  are  de¬ 
scribed.  In  section  3,  results  are  presented  showing  how 
drift  and  variability  are  induced  by  a  change  in  thermal 
damping.  Some  related  issues  are  addressed  in  section 
4.  In  section  5,  major  conclusions  are  presented  and  the 
implications  for  atmosphere -ocean  coupled  modeling 
are  discussed. 

2.  The  model  and  the  thermal  forcing  condition 

This  study  employs  the  Pacanowski  et  al.  (1991) 
version  of  the  Bryan-Cox-Semtner  OGCM,  which  is 
based  on  the  work  of  Bryan  (1969).  The  model  has  a 
width  of  60°  and  a  length  of  144°  extending  from  72°S 
to  72°N.  The  Southern  Hemisphere  includes  a  modeled 
Antarctic  Circumpolar  Current  (ACC)  passage  from 


44°  to  60°S  and  a  sill  2350  m  deep  in  the  model’ s  Drake 
Passage;  elsewhere  it  is  flat  bottomed.  This  sill  is  in¬ 
tended  to  provide  a  realistic  ACC  transport  through  the 
setup  of  a  bottom  pressure  difference  between  the  two 
sides  of  the  passage  (Holland  1973;  Gill  and  Bryan 
1971;  Cai  1994;  Cai  and  Baines  1996).  Otherwise,  the 
flat  bottom  would  lead  to  an  unrealistically  large  wind- 
driven  ACC.  The  horizontal  grid  spacing  is  4°  latitude 
by  4°  longitude.  The  model  has  12  levels  in  the  vertical, 
and  the  depth  distribution  is  listed  in  Table  1.  The 
model  uses  the  Cox  (1987)  parameterization  to  com¬ 
pute  vertical  diffusion  and  convection  implicitly.  The 
enhanced  vertical  diffusivity  in  regions  of  static  insta¬ 
bility  is  set  at  105  cm2  s“\  which  is  the  convective 
adjustment  in  the  model. 

Suppose  a  restoring  spinup  is  carried  out  in  which 
the  model  SST  is  restored  to  Ta,  and  the  steady-state 
solution  has  an  SST  field  Ts  and  a  heat  flux  field  Fspin. 
They  satisfy 

F&pin  =  KH(Ta-Ts).  (1) 

Here  KH  is  the  damping  rate  during  the  spinup.  The 
thermal  forcing  is  then  switched  to 

Qao  =  Kr(Tr-  r.)  (2) 

to  achieve  a  weaker  damping.  Here  Kr  is  the  new  damp¬ 
ing  rate  (Kr  <  KH ),  Tr  is  the  new  restoring  temperature, 
and  Tx  is  the  temperature  of  the  uppermost  level. 

The  Tr  field  is  diagnosed  from  the  spinup  state  by 
setting 

Kr(Tr  -  Ts)  =  Fspin  (3) 

so  that  at  the  moment  of  the  switching  only  the  damp¬ 
ing  rate  changes:  the  heat  flux  is  exactly  the  same  as 
that  prior  to  the  switch.  As  in  a  flux-corrected  coupled 
model,  the  heat  flux  is  subsequently  allowed  to  change. 
The  simplicity  of  such  a  treatment  of  the  forcing  con¬ 
dition  allows  an  examination  of  the  effect  of  a  change 
in  the  rate  of  thermal  damping  alone . 

The  physics  behind  this  treatment  of  surface  thermal 
forcing  can  be  realized  from  a  simple  coupled  system 


Table  1.  Distribution  of  vertical  levels. 


Level 

Thickness  (m) 

Depth  of  T ,  S  (m) 

1 

25.0 

12.5 

2 

25.0 

37.5 

3 

40.0 

70.0 

4 

70.0 

125.0 

5 

110.0 

215.0 

6 

200.0 

370.0 

7 

330.0 

635.0 

8 

450.0 

1025.0 

9 

650.0 

1575.0 

10 

900.0 

2350.0 

11 

900.0 

3250.0 

12 

900.0 

4150.0 
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consisting  of  a  simple  atmosphere  (Schopf  1983;  Sau- 
sen  et  al.  1988;  Zhang  et  al.  1993;  Cai  and  Greatbatch 
1995)  and  a  mixed  layer  ocean.  The  heat  balance  for 
the  atmosphere  is 

0  =  -KH{J.-Tx)-K'rTa  +  Qa  (4) 
and  for  the  mixed  layer  ocean  is 

C0  =  KH(Ta  —  Tx)  +  Q0.  (5) 

at 

Equation  (4)  expresses  a  balance  between  heat  transfer 
from  the  ocean  to  the  atmosphere  given  by  the  -KH(Ta 
-  Ti)  term,  heat  loss  to  space  given  by  -  K'rTa ,  and 
heat  sources  due  to  solar  input  and  the  divergence  of 
the  atmospheric  heat  transport  given  by  Qa.  Here,  K[ 
is  the  atmospheric  radiative  feedback  parameter,  C0  is 
the  specific  heat  capacity  of  the  mixed  layer  ocean,  and 
Q0  is  the  oceanic  heat  source. 

Defining  a  temperature 

Tr  =  Qa/K'r,  (6) 


we  can  write  (4)  as 


0  =  -K„(Ta  -  TO  +  Kfr(Tr  -  Ta).  (7) 

From  (7),  expressing  Ta  in  terms  of  other  quantities 
and  substituting  into  ( 5  ) ,  we  have 

c0^  =  Kr{Tr-  ro  +  a.  (8) 

at 


where 


Kr  = 


KhK’ 

(Kh  +  K'r)  ‘ 


(9) 


A  comparison  of  (5)  and  (8)  indicates  that  (2)  represents 
a  simple  coupling  [see  Zhang  et  al.  ( 1993 )  for  a  detailed 
discussion] .  Since  A''  is  a  factor  of  about  20  less  than  KH 
(Dickinson  1981),  it  follows  from  (9)  that  Kr  is  also 
much  smaller  than  KH.  Thus,  switching  from  restoring 
forcing  to  (2)  leads  to  a  weaker  thermal  damping.  In  (2) 
Tr  is  kept  constant.  This  means  that  Qa  is  fixed  [see  Eq. 
(6)].  The  physical  implication  is  that  the  atmospheric 
heat  transport  under  the  new  forcing  is  exactly  the  same 
as  that  implied  by  the  steady  state  of  the  spinup. 

Defining  a  ratio 


a  =  Kri Kh  (a  <  1 ),  (10) 

which  measures  the  change  in  thermal  damping,  we  can 
rewrite  the  new  forcing  [i.e.,  (2)],  using  (1)  and  (3), 
as 


Kr  ( Tr  T\  ) 


=  Kr(Tr  -  Ts)  -  Kr(Ta  -  Ts)  +  Kr(Ta  -  Tx) 

=  (1  -  a)Fspin  +  aK„(Ta  -  Tx).  (11) 


Thus,  the  new  forcing  is  equivalent  to  a  part  of  the 
constant  diagnosed  heat  flux  plus  a  weaker  restoring 
forcing. 

3.  Model  results 

a.  Change  of  thermal  damping  from  a  spinup  under 

restoring  conditions  on  both  temperature  and 

salinity 

The  OGCM  is  spun  up  by  restoring  the  uppermost 
level  temperature  and  salinity  to  a  zonally  uniform  tem¬ 
perature  and  salinity,  given  by  the  profiles  shown  in 
Fig.  1.  These  approximate  the  zonally  averaged  value 
of  the  observed  fields.  The  e-folding  restoring  time  is 
30  days  for  salinity  and  15  days  for  temperature,  the 
latter  corresponding  to  a  thermal  damping  rate  KH  of 
79  W  m-2  °C~1 .  We  choose  a  longer  timescale  for  sa¬ 
linity  to  produce  a  freshwater  flux  with  more  realistic 
magnitudes.  This  prevents  a  collapse  of  the  overturning 
circulation  upon  switching  to  mixed  boundary  condi¬ 
tions  (Tziperman  et  al.  1994;  Cai  1996a).  The  model 
is  also  subject  to  the  wind  stress  used  by  Bryan  ( 1987) . 
Only  during  the  first  1 1  000  years  (88  000  years  at  the 
bottom  level)  of  the  spinup  is  the  acceleration  of  Bryan 
( 1984)  used.  By  then,  the  spinup  reaches  a  statistically 
steady  state.  Then,  the  integration  is  extended  for  an¬ 
other  1000  years  without  the  acceleration  for  a  final 
steady  state.  The  overturning  circulation  is  shown  in 
Fig.  2  and  is  qualitatively  similar  to  that  in  the  real 
Atlantic.  It  features  a  northern  sinking  cell,  similar  to 
the  North  Atlantic  Deep  Water  formation  (NADWF) 
cell  and  an  Antarctic  Circumpolar  Current  convective 
cell.  The  model  ACC  during  the  spinup  reaches  102  Sv 
(Sv  =  106  m3  s'"1). 

In  run  1  (Table  2),  the  thermal  damping  rate  is 
switched  to  Kr  —  19.75  W  m“2  °C_1  (corresponding  to 
a  damping  timescale  of  60  days),  with  the  forcing 
given  by  (2)  and  Tr  diagnosed  from  the  steady  state  in 
terms  of  (3) .  In  this  case,  a  =  0.25,  and  from  ( 1 1 )  the 
new  forcing  consists  of  75%  of  the  diagnosed  heat  flux 
and  25%  of  the  restorative  forcing  used  in  the  spinup. 
The  restoration  for  salinity  is  maintained.  It  is  relevant 
to  note  that  a  flux  correction  for  coupled  models  is  usu¬ 
ally  diagnosed  by  averaging  the  flux  fields  over  a  pe¬ 
riod  of  several  years.  This  does  not  ensure  an  identical 
heat  flux  at  the  moment  of  the  switch  and  may  exert 
effect  on  the  subsequent  model  solution.  In  order  to 
eliminate  this  effect  in  the  present  study,  such  an  av¬ 
erage  is  not  carried  out. 

The  model  is  integrated  for  a  further  1 1  000  years. 
Strong  drift  occurs  as  can  be  seen  from  Fig.  3,  which 
shows  the  time  series  of  the  overturning  sampled  at 
latitude  60°N,  depth  2350  m.  Unless  otherwise  stated, 
single-point  time  series  of  overturning  are  all  taken  at 
this  location.  Several  features  emerge.  First,  the  new 
steady-state  drifts  from  the  spinup  state,  and  the  larg¬ 
est  timescale  for  the  drift  reaches  several  thousands  of 
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a),  Relaxation  Temperature 


b),  Relaxation  Salinity 


Fig.  1 .  (a)  Surface  relaxation  temperature  (in  °C)  and  (b)  surface 
relaxation  salinity  (in  parts  per  thousand). 


years.  Second,  some  variabilities,  apparently  ther¬ 
mally  driven,  are  induced  at  initial  stage.  When  the 
thermal  damping  is  weakened  further,  to  Kr  =  9.875 
W  mT2  °C_1  ( a  =  0.125),  giving  rise  to  run  2  (Table 
1),  regular  thermal  oscillations  with  a  period  of  26 
years  occur  (Fig.  4).  The  characteristics  of  the  oscil¬ 
lation  are  given  in  Greatbatch  and  Zhang  (1995)  and 
Cai  et  al.  (1995)  and  will  not  be  repeated  here.  As  in 
these  studies,  the  oscillation  reproduced  some  impor¬ 
tant  features  of  the  variability  seen  in  the  Geophysical 


Fluid  Dynamical  Laboratory  (GFDL)  coupled  model 
(Del worth  et  al.  1993). 

An  examination  reveals  that  the  location  of  the  max¬ 
imum  overturning  of  run  1  shifts  away  from  that  of  the 
spinup  and  the  maximum  decreases  by  about  0.5  Sv. 
This  drift  significantly  modifies  the  temperature  field, 
especially  the  SSTs  (Fig.  5a).  This  modification  is  as¬ 
sociated  with  the  Tr  field  (Fig.  5b)  and  changes  in  the 
convective  pattern  (Figs.  5c  and  5d).  By  convection 
pattern  we  mean  the  location  and  intensity  (the  pene¬ 
tration  depth)  of  the  convection.  In  order  to  maintain 
the  convection  pattern  of  the  spinup,  Tr  over  the  con¬ 
vective  regions  is  “corrected”  to  be  fairly  cold  in  order 
to  allow  heat  loss.  Once  the  convection  pattern  (loca¬ 
tion  and  intensity)  changes,  this  cold  bias  becomes  the 
cause  for  drift.  For  example,  under  the  weaker  damp¬ 
ing,  convection  in  the  northeastern  comer  weakens,  and 
some  even  cease,  but  the  low  value  of  Tr  in  that  region 
still  allows  large  heat  loss.  This  leads  to  a  drop  in  SST. 
A  similar  process  takes  place  in  the  eastern  side  of  the 
model  Weddell  Sea,  and  the  cold  water  appears  to 
spread  along  the  model  ACC.  By  contrast,  convection 
off  the  northern  western  boundary  intensifies,  where  Tr 
is  not  low  enough  to  remove  the  associated  heat  trans¬ 
ported  from  below  by  the  intensified  convection.  This 
leads  to  a  strong  warming  there.  This  process  is  also 
present  in  run  2,  where  the  convection  pattern  changes 
periodically  (Fig.  6). 

Instability  of  the  convection  pattern  of  the  spinup 
state  is  a  major  feature  associated  with  the  drift.  This 
is  supported  by  Fig.  5a,  which  shows  that  the  largest 
drifts  occur  over  convective  regions.  What  causes  the 
instability?  In  a  study  on  the  effect  of  a  change  in  saline 
damping  using  a  low -resolution  model  with  realistic 
topography,  Moore  and  Reason  (1993)  discussed  the 
source  of  drift  and  variability  in  their  model.  The  so- 
called  equilibrium  solution  of  the  spinup  via  relaxation 
to  a  given  climatology  is  not  a  perfect  steady-state  so¬ 
lution  of  the  model  equations.  In  areas  characterized  by 
sharp  oceanic  fronts  and  high  convective  activity,  the 
model,  due  to  the  limitation  of  model  dynamics,  cannot 
maintain  the  climatology.  Therefore,  in  these  areas  the 


Meridional  Overturning  (Sv)  From  Spin-up 


Fig,  2.  Overturning  streamfunction  (in  Sverdrups)  from  the  spinup 
(run  RR)  under  restoring  boundary  conditions  for  temperature  and 
salinity. 
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Table  2.  Details  of  main  model  experiments.  The  thermal  forcing  condition  is  either  a  restoring  with  the  restoring  time  (in  days)  indicated 
in  parentheses  or  that  described  by  Eq.  (11),  that  is,  a  restoring  plus  a  portion  of  a  diagnosed  flux.  The  saline  forcing  condition  is  either  a 
restoring  or  diagnosed  flux.  In  the  thermal-only  cases,  the  salinity  forcing  is  zero. 


Run 

Initial 

condition 

Thermal  forcing 

Saline  forcing 

Characteristics 

RR 

Rest 

Restoring  (15  d) 

Restoring  (30  d) 

Steady 

1 

RR 

Restoring  (60  d)  +  Flux 

Restoring  (30  d) 

Steady,  drift  from  RR 

2 

RR 

Restoring  (120  d)  +  Flux 

Restoring  (30  d) 

Oscillatory,  drift  from  RR 

RF 

RR 

Restoring  (15  d) 

Flux 

Steady 

3 

RF 

Restoring  (120  d)  +  Flux 

Flux 

Oscillatory,  drift  from  RF 

R0 

Rest 

Restoring  (15  d) 

(Thermal  only) 

Steady 

4 

Rest 

Restoring  (30  d)  4-  Flux 

(Thermal  only) 

Steady,  drift  from  R0 

5 

Rest 

Restoring  (45  d)  +  Flux 

(Thermal  only) 

Oscillatory,  drift  from  R0 

model  is  systematically  biased  away  from  the  clima¬ 
tology,  and  the  surface  heat  and  salinity  fluxes  tend  to 
oscillate  about  some  mean  value.  This  idea  was  sup¬ 
ported  by  time  series  of  the  global  average  fluxes, 
which  show  that  as  the  imperfect  steady  state  is  ap¬ 
proached,  the  global  average  oceanic  surface  heat  and 
salt  fluxes  possess  a  significant  oscillatory  component. 
Similar  time  series  from  our  spinup  confirm  theirs.  The 
oscillatory  feature  continues  to  be  present  even  if  our 
model  is  integrated  for  a  further  10  000  years. 

During  the  spinup,  SSTs  are  forcefully  4  ‘nudged’  ’  to 
the  prescribed  values,  in  a  similar  way  to  some  data 
assimilation  studies.  In  areas  where  the  model  is  not 
compatible  with  the  climatology  it  tries  to  resist  the 
restoring  (nudging)  conditions  at  the  surface  through 
the  action  of  dynamic  adjustment  [see  Moore  and  Rea¬ 
son  (1993 )  and  the  references  therein] .  The  adjustment 
process,  however,  never  reaches  completion  because 
the  restoring  conditions  constantly  damp  it.  The  restor¬ 
ing  boundary  conditions  therefore  act  as  a  “brake”  on 
the  dynamic  adjustment  process  and  suppress  the  “in¬ 
compatibility.”  This  leads  to  intermittent  convections 
of  the  spinup  steady  state. 

Once  the  surface  restoring  condition  on  temperature 
is  switched  to  a  weak  restoring  condition  ( weak  damp¬ 
ing),  temperature  anomalies  in  areas  with  incompati¬ 
bility  are  freer  to  evolve.  In  other  words,  the  internal 
model  dynamics  is  allowed  to  play  a  greater  role  in 
determining  the  convection.  The  intermittent  convec¬ 
tions  of  the  spinup  state  change  their  location  and  pen¬ 
etration  depth.  This  leads  to  the  drift. 

If  the  drift  is  substantial,  it  influences  the  convection 
pattern  in  turn  and  produces  oscillations,  a  periodic  pro¬ 
cess  in  which  the  anomaly  associated  with  the  drift  and 
the  change  in  the  convection  pattern  each  become  the 
cause  for  the  other.  This  can  be  understood  via  the  ad- 
vective  process  (Huang  and  Chou  1994;  Greatbatch 
and  Zhang  1995;  Cai  et  al.  1995).  Suppose  that  at  a 
location  a  negative  SST  anomaly  develops  due  to  drift. 
The  anomaly  enhances  convection,  and  more  subsur¬ 
face  heat  is  transmitted  to  the  surface.  If  the  surface 
cooling  is  not  strong  enough  to  remove  the  additional 
heat,  a  positive  SST  anomaly  is  produced.  The  positive 


SST  anomaly  then  reduces  the  convection,  and  SST 
drops,  owing  not  only  to  the  weakened  subsurface  heat¬ 
ing  but  also  to  the  surface  cooling.  Consequently,  a 
negative  SST  anomaly  develops.  In  this  way,  the  cycle 
repeats  itself,  provided  that  the  damping  is  weak.  Al¬ 
ternatively,  the  oscillation  can  be  understood  via  the 
thermohaline  adjustment  considered  by  Wajsowicz  and 
Gill  ( 1986)  in  terms  of  boundary  waves.  This  process 
was  recently  revisited  by  Winton  (1996)  and  Great- 
batch  and  Paterson  ( 1996).  These  authors  find  that  an 
anomalous  northward  flow  associated  with  a  stronger 
overturning  pumps  a  heat  anomaly  into  the  northeastern 
Conner.  The  anomaly  then  propagates  westward  along 
the  northern  boundary  by  means  of  a  boundary  wave. 
Much  of  the  time  in  an  oscillation  cycle  is  spent  at  the 
northern  boundary,  where  the  weak  vertical  stratifica¬ 
tion  slows  the  propagation  of  boundary  waves,  leading 
to  the  emergence  of  an  interdecadal  timescale. 

In  assessing  the  relevance  of  our  model  results  to 
oscillations  in  a  coupled  model,  we  note  that  in  a  cou¬ 
pled  model  atmospheric  changes  also  take  place.  These 
atmospheric  changes  are  likely  to  play  a  role  in  gen¬ 
erating  oscillations.  This  is  proposed  by  Griffies  and 
Tziperman  (1995)  in  a  box  model  study  as  a  cause  for 
the  oscillation  in  the  GFDL  coupled  model  (Del worth 
et  al.  1993).  Their  result  suggests  that  the  variability 
in  the  coupled  model  is  externally  driven  rather  than 
being  an  internal  oceanic  process,  as  discussed  in  the 
last  paragraph.  Clearly,  there  is  a  need  to  reconcile 
these  results.  It  is  likely  that  in  a  coupled  system  these 
processes  coexist.  This  is  an  important  aspect  that 
needs  further  study. 

Returning  to  the  diift  and  the  oscillation  in  the  pres¬ 
ent  study,  we  note  that  our  results  are  in  sharp  contrast 
to  those  of  Sausen  et  al.  (1988),  who  reported  that 
when  a  heat  flux  correction  was  employed  and  the  res¬ 
toration  for  salinity  was  maintained,  no  drift  occurred. 
It  is  important  to  point  out  that  their  model  spinup  was 
in  a  perfect  steady  state;  hence,  there  was  no  internal 
variability  that  could  grow  upon  the  change  in  the  ther¬ 
mal  damping.  In  our  model  drifts  occur  regardless  of 
the  forcing  condition  for  salinity.  This  can  be  seen  from 
the  following  manipulation: 
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a),  Time  Series  of  Overturning,  Run  1,  Years  0  -  1000 


Time  (years) 


b).  Time  Series  of  Overturning,  Run  1,  Years  1000  -  4000 


Fig.  3.  Time  evolution  of  the  northern  overturning  (in  Sverdrups) 
at  a  location  of  64°N,  2350  m  for  run  1,  that  is,  with  a  thermal  damp¬ 
ing  of  19.75  W  m"2  K_I  and  a  restoration  on  surface  salinity,  for  (a) 
the  first  1000  years  and  (b)  the  period  between  year  1000  and  year 
5000. 


Qao  =  Kr(Tr  -  Ti)  =  Kr(Tr  -  Ts)  +  Kr(Ts  -  T{) 

-  Fspin  +  Kr(Ts  —  7\).  (12) 

Equation  (12)  states  that  although  the  heat  flux  is  the 
same  as  in  the  spinup,  the  ocean  will  drift  from  Ts. 
Further,  as  Kr  decreases,  the  SST  is  freer  to  evolve  from 
the  modeled  climatology.  An  examination  reveals  that 
in  both  run  1  and  run  2,  salinity  (especially  SSS)  also 
drifts  from  that  of  the  spinup,  in  response  to  changes 
in  convection  pattern. 


b.  Change  of  thermal  damping  from  a  spinup  under 

mixed  boundary  conditions 

In  order  to  further  test  this  hypothesis  of  the  in¬ 
evitability  of  the  drift  we  have  also  carried  out  an¬ 
other  experiment,  designated  run  3.  It  starts  from  a 
steady  state  of  a  second  spinup  RF  (Table  2),  which 
is  forced  by  mixed  boundary  conditions,  with  a  salt 
flux  diagnosed  also  from  RR.  The  model  under  mixed 
boundary  conditions  has  similar  overturning  circu¬ 
lations  to  those  shown  in  Fig.  2  owing  to  the  weak 


b),  Time  Series  of  Overturning,  Run  2.  Years  4000  -  5000 


Fig.  4.  Time  evolution  of  the  northern  overturning  (in  Sverdrups) 
at  64°N,  2350  m  for  run  2,  that  is,  with  a  thermal  damping  of  9.875 
W  m-2  KT1  and  a  restoration  on  surface  salinity,  for  (a)  period  be¬ 
tween  year  0  and  year  5000  and  (b)  a  magfication  covering  the  period 
between  year  4000  and  year  5000. 
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a),  SST  Drift,  Run  1 


0  30E  60E 

c)t  Convection  Pattern,  Run  RR 


b),  Tr  Fietd,  Run  1 
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d)r  Convection  Pattern,  Run  1 


Fig.  5.  (a)  SST  drift  for  run  1 .  The  drift  is  the  difference  between  an  average  over 
the  last  10  years  of  run  1  and  the  steady  state  of  run  RR.  (b)  Tr  field  used  to  drive 
run  1.  (c)  and  (d)  convection  patterns  (i.e.,  the  locations  and  intensity  of  convections) 
of  the  steady  state  of  run  RR  and  run  1,  respectively.  The  intensity  is  shown  by  the 
penetration  depth  in  terms  of  model  levels,  indicated  by  the  legend. 


restoration  for  SSS  (Tziperman  et  al.  1994;  Cai 
1996a).  The  thermal  forcing  condition  is  then 
switched  to  one  the  same  as  that  of  run  2,  and  the 
flux  forcing  for  salinity  remains.  Drift  and  variability 
occur  (Fig.  7)  similar  to  those  in  run  2  but  with  a 
period  of  27  years.  The  oscillation  amplitude  de¬ 


creases  somewhat  owing  to  the  ‘‘braking”  effect  of 
salinity  because  SSS  is  free  to  respond.  Figure  8 
shows  the  difference  fields  of  SST  and  SSS  between 
the  mean  over  the  last  10  steady  oscillation  cycles 
and  the  initial  state.  Indeed,  substantial  SSS  drift  oc¬ 
curs,  even  though  only  the  thermal  damping  alters. 
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a).  Peak  b),  Mid.  Peak-to-Trough 
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Fig.  6.  Convection  pattern  at  a  quarter  period  interval  for  run  2  corresponding  to 
times  when  the  overturning  in  Fig.  3b  is  at  (a)  a  maximum  and  (b)  at  a  minimum. 


c.  Change  of  thermal  damping  from  a  spinup  under  a 
thermal  forcing  only 

We  have  also  carried  out  experiments  in  which  sa¬ 
linity  is  kept  constant  and  uniform  at  34  ppt.  Initially, 
the  model  is  spun  up  to  a  steady  state  by  restoring 
the  model  surface  temperature  to  that  given  in  Fig. 
la  using  a  KH  =  79  W  m"2  °C_1 .  This  run  is  referred 
to  as  run  R0  (Table  2).  Runs  4  and  5  start  from  the 


steady  state  of  R0,  with  a  Kr  of  39.5  (a  =  0.5)  and 
26.3  (a  =  0.33)  W  m~2  °C~1,  corresponding  to  a 
damping  timescale  of  30  and  45  days,  respectively 
(Table  2).  The  model  is  integrated  for  a  further 
11  000  years.  Drift  takes  place  in  both  runs,  and  in 
run  5  interdecadal  oscillations  with  a  period  of  27.1 
years  persist  (Fig.  9).  Comparing  runs  3  and  5,  we 
note  that  persistent  oscillations  are  more  easily  gen¬ 
erated  in  run  5.  This  is  because  in  run  5  the  braking 
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a).  Time  Series  of  Overturning,  Run  3 
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Fig.  7.  Time  evolution  of  the  northern  overturning  (in  Sverdrups) 
at  64°N,  2350  m  for  run  3,  which  starts  from  a  steady  state  of  a  run 
(run  RF)  under  mixed  boundary  conditions.  From  the  steady  state  of 
RF,  the  thermal  forcing  condition  is  diagnosed  using  a  thermal  damp¬ 
ing  rate  of  9.875  W  m~2  K-t:  (a)  the  period  between  year  0  and  year 
5000  and  (b)  a  magnification  covering  the  period  between  year  4000 
and  year  5000. 


effect  of  salinity  is  absent  [see  Greatbatch  and  Zhang 
(1995)  and  Cai  et  al.  (1996b)  for  a  discussion]. 

4.  Discussion 

Before  we  proceed  to  discuss  the  relevance  of  these 
model  results  to  the  coupled  atmosphere -ocean  mod¬ 
eling,  it  is  useful  to  address  several  issues.  First,  during 
our  spinups  (runs  RR,  RF,  or  R0)  we  use  a  thermal 


damping  rate  of  79  W  m  2  °C  -1 ,  a  rate  about  twice  that 
commonly  used.  To  study  the  possible  effect  of  this 
large  initial  damping,  we  have  conducted  runs  that  con¬ 
tinue  from  spinups  under  a  thermal  damping  rate  of 
39.5  W  m~2  °C  _1 .  We  find  that  although  the  magnitude 
of  the  drift  varies  with  the  extent  of  the  change  in  ther¬ 
mal  damping,  the  major  features  described  above  ap¬ 
pear.  That  is,  as  the  thermal  damping  weakens,  drift 
occurs;  as  the  thermal  damping  weakens  further,  inter¬ 
nal  variability  takes  place. 

Second,  the  design  of  the  experiments  may  favor  a 
substantial  drift  and  a  vigorous  oscillation  when  the 
thermal  damping  weakens.  In  our  study,  a  simple 
model  is  run  with  zonally  averaged  restoring  thermo¬ 
haline  conditions,  which  are  derived  from  observations. 
Such  a  simple  approach  renders  that  one  has  no  way  of 
knowing  a  priori  whether  these  restoring  conditions 
could  be  compatible  with  the  OGCM.  It  is  reasonable 
to  argue  that  models  with  realistic  basin  configuration 
and  bottom  topography  would  have  a  better  likelihood 
of  being  compatible  with  the  restoring  conditions  since 
these  restoring  conditions  are  based  on  a  realistic  sys¬ 
tem.  Although  this  may  be  true,  we  believe  that  incom¬ 
patibility  generally  exists.  We  know  that  the  observed 
forcing  conditions  are  not  perfect;  even  if  the  OGCM 
is  state  of  the  art,  the  model  physics  do  not  fully  rep¬ 
resent  that  in  the  real  world.  Deficiencies  in  model  dy¬ 
namic  processes,  either  absent,  poorly  parameterized, 
or  not  accommodated  by  the  model  resolution,  mean 
that  there  is  incompatibility.  Therefore,  drift  and  oscil¬ 
lation  will  occur  once  the  incompatibility  is  allowed  to 
manifest  itself,  although  the  intensity  of  these  features 
varies  with  the  extent  of  the  incompatibility. 

Third,  in  the  presence  of  any  incompatibility,  con¬ 
vection  in  the  spinup  state  is  intermittent.  That  the  in¬ 
stability  of  a  convection  pattern  plays  a  significant  role 
in  the  generation  of  drift  and  oscillations  raises  a  ques¬ 
tion  as  to  how  the  model  results  depend  upon  the  model 
convection  scheme.  To  address  this  question,  we  have 
repeated  run  RR  (a  spinup)  and  run  1  using  the  con¬ 
vection  scheme  of  Rahmstorf  (1993).  The  advantage 
of  the  scheme  over  that  of  Cox  (1987)  is  that  it 
achieves  complete  mixing  when  a  vertical  instability 
occurs.  We  found  that  upon  switching  to  weaker  ther¬ 
mal  damping,  the  model  behaves  in  a  way  similar  to 
that  in  run  1 . 

Finally,  as  commented  in  section  3a,  a  flux  correc¬ 
tion  for  a  coupled  model  is  usually  diagnosed  by  av¬ 
eraging  the  flux  fields  over  a  period  of  several  years. 
In  order  to  ensure  an  identical  heat  flux  at  the  moment 
of  the  switch,  this  is  not  carried  out  in  the  experiments 
described  above.  To  test  the  possible  effect,  we  have 
repeated  run  1  using  a  heat  flux  diagnosed  by  averaging 
over  the  final  10  years  of  the  spinup.  The  general  time- 
dependent  behavior  and  final  state  of  the  model  hardly 
differ  from  those  in  the  original  run  1 .  Thus,  we  con¬ 
clude  that  the  effect  of  averaging  the  heat  flux  is  neg- 
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Fig.  8.  (a)  SST  and  (b)  SSS  drift  of  a  mean  state  of  run  3  from  the  steady  state 
of  run  RF.  The  mean  state  is  averaged  over  the  last  10  steady  oscillation  cycles, 
that  is,  over  the  period  between  year  10  730  and  year  1 1  000. 


ligible,  and  it  is  the  change  in  thermal  damping  that  has 
predominant  effects. 

5.  Conclusions  and  implications  for  coupled 

modeling 

In  atmosphere -ocean  coupling,  a  precoupling 
spinup  of  an  OGCM  via  Haney  restoration  to  the  pres¬ 
ent  day  climatology  is  commonly  carried  out  to  provide 
initial  conditions  and,  in  many  cases,  to  construct  flux 
correction  fields  for  coupled  models.  Haney  restoration 
represents  a  strong  damping,  and  it  is  now  proven  be¬ 
yond  doubt  that  upon  switching  to  an  interactive  at¬ 
mosphere  the  thermal  damping  changes  from  strong, 
constant,  and  uniform  damping  to  one  that  is  generally 
weaker  and  time  and  scale  dependent.  In  an  intercom¬ 
parison  of  the  surface  energy  budget  for  a  perpetual 
July  among  19  atmospheric  general  circulation  models 
( AGCMs),  Randall  et  al.  ( 1992)  find  that  the  average, 
over  the  global  domain  (including  the  land  area)  and 
over  all  the  19  AGCMs,  rate  of  thermal  damping  in 
response  to  a  global-scale  SST  anomaly  is  about  1.9 
W  m'2  K_1 .  When  the  land  area  is  excluded,  the  rate 
may  be  larger  but  it  is  still  much  smaller  than  that  dur¬ 
ing  the  precoupling  spinup.  For  example,  the  CSIR09 
AGCM,  when  forced  in  the  same  way  as  the  other  19 
AGCMs,  gives  a  value  of  2.1  W  m-2  K  1  over  a  global 
domain;  when  the  land  area  is  excluded,  the  value  is 
3.4  W  m~2  K-1,  much  smaller  than  that  in  the  spinup. 
These  AGCM  results  reinforce  the  need  for  studying 


the  effect  of  a  change  in  thermal  damping  when  cou¬ 
pling  occurs. 

The  model  results  of  the  present  study  suggest  that 
climate  drift  will  occur  in  either  a  flux-corrected  or  an 
uncorrected  coupled  model.  The  weakening  in  the  ther¬ 
mal  damping  allows  the  internal  model  dynamics  to 
play  a  greater  role  in  determining  the  location  and  in¬ 
tensity  of  convections,  which  are  generally  very  dif¬ 
ferent  from  those  during  the  spinup.  This  difference 
causes  a  drift,  a  feature  also  seen  in  flux-corrected  cou¬ 
pled  models  (e.g.,  GOF). 

The  timescale  of  the  drift  suggests  that  climate  sce¬ 
nario  experiments  should  be  initiated  long  after  switch¬ 
ing  to  the  coupled  mode,  that  is,  after  the  ocean  has 
settled  to  a  new  steady  state,  so  that  the  scenario  event 
is  not  contaminated  by  drift.  This  point  is  important 
only  if  the  magnitude  of  the  drift  is  of  comparable  order 
to  the  signal  being  detected. 

As  the  thermal  damping  weakens,  regular  thermal 
oscillations  take  place.  Oscillations  at  interdecadal 
timescales  are  produced  in  both  the  GFDL  (Delworth 
et  al.  1993)  and  CSIRO  (GOF)  coupled  model.  Great- 
batch  and  Zhang  (1995)  and  Cai  et  al.  (1995)  have 
shown  that  the  oscillation  in  the  GFDL  coupled  model 
resembles  well  the  thermally  driven  oscillations  in  their 
models.  This  suggests  that  oscillations  in  coupled  mod¬ 
els  are  likely  to  be  induced  by  the  weakening  in  thermal 
damping  associated  with  the  coupling,  although  we 
note  that  Griffies  and  Tziperman  ( 1995 )  have  proposed 
a  new  explanation  for  the  generation  of  the  oscillation 
described  by  Delworth  et  al.  (1993). 
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a),  Time  Series  of  Northern  Overturning,  Run  4 


b),  Time  Series  of  Northern  Overturning,  Run  5 


Time  {years) 

Fig.  9.  Time  evolution  of  the  northern  overturning  (in  Sverdrups) 
at  64°N,  2350  m  for  runs  4  and  5.  Both  start  from  the  steady  state  of 
a  run  (run  R0)  under  a  thermal  forcing  only.  The  thermal  forcing 
conditions  are  diagnosed  from  the  steady  state  using  damping  rate 
39.5  (run  4)  and  26.33  (run  5)  W  m-2  K~‘. 


Most  of  the  precoupling  spinups  for  atmosphere - 
ocean  coupled  models  are  carried  out  by  restoring 
model  SST  and  SSS  to  Levitus  climatologies.  These 
spinups  will  have  problems  obtaining  a  perfect  steady 
state,  especially  when  the  model  resolution  is  coarse. 
In  these  models,  it  is  not  clear  whether  the  given  cli¬ 
matologies  can  be  maintained  by  model  dynamics.  As 
pointed  out  by  Moore  and  Reason  (1993),  in  areas 
characterized  by  sharp  fronts  in  the  climatologies,  mod¬ 
els  experience  problems  maintaining  strong  fronts  such 
as  those  observed  because  the  large  horizontal  eddy 


diffusivities  in  the  model  constantly  try  to  eliminate 
them.  This  means  that  in  these  regions,  the  climatolo¬ 
gies  are  not  compatible  with  the  model.  In  the  presence 
of  such  incompatibility,  the  model  convection  in  the 
spinup  steady  state  is  intermittent,  and  the  atmospheric 
heat  and/or  freshwater  transport  implied  in  the  diag¬ 
nosed  heat  and  freshwater  fluxes  are  not  a  perfect  match 
with  the  oceanic  heat /freshwater  transports.  Cai 
( 1 995 )  has  demonstrated  that  such  a  mismatch  can  gen¬ 
erate  interdecadal  oscillations  when  the  diagnosed 
fluxes  are  used  to  force  the  model. 

A  flux  correction  diagnosed  from  such  a  spinup  bears 
the  imprint  of  the  incompatibility,  which  is,  upon  cou¬ 
pling,  permanently  present  in  the  subsequent  simula¬ 
tions.  As  discussed  above,  this  may  become  the  cause 
of  drift  and  variability.  On  the  other  hand,  there  is  no 
dynamical  reason  why  in  coupled  models  the  heat  flux 
from  the  model  atmosphere  to  the  model  ocean  should 
always  match  the  oceanic  dynamics.  For  example,  the 
heat  flux  is  influenced  by  many  atmospheric  processes, 
which  are  weakly  dependent  or  even  completely  inde¬ 
pendent  of  the  ocean  (Cai  et  al.  1995;  Cai  1995).  At 
this  stage,  however,  it  is  not  entirely  clear  what  gen¬ 
erates  the  oscillations  in  the  existing  fully  coupled  mod¬ 
els  (e.g.,  Delworth  et  al.  1994;  GOF).  In  our  model, 
this  incompatibility  is  carried  over  to  the  new  system 
via  the  constant  flux  term  in  (11)  (see  §2),  and  the 
incompatibility  then  drives  the  drift  and  the  oscillation. 

That  the  convection  pattern  changes  with  the  thermal 
damping  raises  questions  regarding  the  validity  of  the 
flux  correction  as  diagnosed  from  a  restoring  spinup.  It 
is  tempting  to  argue  that  the  correction  should  perhaps 
be  derived  from  a  more  realistic  simulation.  For  ex¬ 
ample,  run  3  is  performed  using  a  more  realistic  damp¬ 
ing  rate,  and  therefore  the  flux  correction  should  be 
diagnosed  from  such  a  run.  This  new  flux  correction 
will  be  rather  different  from  that  diagnosed  from  the 
spinup,  because  it  has  to  maintain  a  different  convec¬ 
tion  pattern.  This  may  reduce  the  drift  somewhat,  but 
since  the  thermal  forcing  condition  of  this  run  is  di¬ 
agnosed  from  the  spinup,  it  still  inherits  the  incompat¬ 
ibility.  Further,  if  the  diagnosed  correction  has  to  adjust 
for  the  maintenance  of  incompatibility-driven  oscilla¬ 
tions,  it  would  need  to  be  time  dependent,  varying  on 
the  timescale  of  the  oscillation.  This  is  not  only  im¬ 
practical  but  also  has  the  effect  of  preconditioning  an 
oscillatory  signal.  It  follows  that  an  ultimate  reduction 
of  incompatibility-induced  drift  and  variability  can 
only  be  achieved  by  improving  model  dynamics  and 
the  quality  of  climatologies. 
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